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We present the results from full QCD simulations with four flavors of light staggered dynamical quarks on 
QCDSP supercomputer. Previous results are reproduced and the simulation reported here yields new results 
consistent with our previous runs. The hadron spectrum obtained with Wilson valence fermions reported here 
will allow us to determine if our earlier conclusions are independent of lattice formalism. 
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1. INTRODUCTION 

It is well known that the number of light quark 
flavors Nf governs the order of the QCD phase 
transition at finite temperature. On the other 
hand the effects of dynamical quarks at zero tem- 
perature have generally been found to be quite 
small when zero-flavor (quenched) and two-flavor 
QCD simulations are compared. One important 
question is whether this will continue as Nf in- 
creases. In the past few years the Columbia group 
has done extensive comparison studies of zero- 
flavor, two-flavor and four-flavor QCD at zero 
temperature with the lattice spacing and vol- 
ume fixed in physical units determined by the p 
mass. What we found was that for four flavors the 
hadron spectrum and chiral condensate were sig- 
nificantly different from zero and two flavors: the 
hadronic spectrum appeared to be parity doubled 
and the chiral condensate was almost four times 
smaller than that for two flavors|l|[^). 

Those calculations were done on the Columbia 
16 Gflop machine, Fermi256. With the new gen- 
eration of Columbia supercomputer, QCDSP, we 
have been able to continue our flavor-dependence 
study. As a first step we need to check the relia- 
bility of our new hardware and software. We have 
reproduced results from several other groups S on 
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small lattices, and we have also tried to reproduce 
our own earlier results. For our previous four- 
flavor calculations we did two simulations with 
dynamical quark mass 0.01 and 0.02. It is natural 
to have at least one more data point to be able to 
do a linear fit to the dynamical results. Another 
interesting question is whether the partial chiral 
symmetry restoration that we observed is just an 
artifact of the lattice formalism. To investigate 
this problem we have also measured the Wilson 
valence spectrum on these lattices generated with 
staggered dynamical quarks. 

2. SIMULATIONS AND RESULTS 

Table [l] lists the run parameters for our four- 
flavor simulations. Run I and III are our previous 
simulations done on Fermi256Q| while II and IV 
are simulations performed on QCDSP. Run IV is 
a check of results from QCDSP against previous 
results from Fermi256. Run II is a new simula- 
tion, adding one more dynamical data point so 
that linear fits can be performed. The four-flavor 
simulations used an exact hybrid Monte Carlo al- 
gorithm, employing the <J> hybrid molecular dy- 
namics algorithm of g, with a Monte Carlo ac- 
cept/reject step. 

Comparing the hadron masses in run III and 
IV in Table [j], we find that the hadron spectrum 
agrees very well statistically between the two ma- 
chines. The expectation values of average pla- 
quette and chiral condensate are also found to 
be consistent. However there is one discrepancy 
between the two calculations which we don't un- 
derstand. In order to achieve comparable accep- 
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Table 1 

Simulation parameters for runs with Nf — 4 and (3 — 5.4 on 16 3 x 32 lattice and corresponding results 
of staggered hadron masses with m va i = md yn - We have chosen the hadron mass fitting range typically 
from 9 to 16 for the it and from 6 to 16 for the other particles. 



N f = 4, = 5.4 


run 


I 


II 


III 


IV 


m dyn a 


0.01 


0.015 


0.02 


0.02 


run length 


4450 


4390 


2725 


3980 


step size 


0.0078125 


0.005 


0.01 


0.005 


steps per trajectory 


64 


100 


50 


100 


acceptance 


0.95 


0.87 


0.99 


0.91 


machine 


Fermi256 


QCDSP 


Fcrmi256 


QCDSP 


7T 


0.292(5) 


0.320(4) 


0.357(3) 


0.356(2) 


P 


0.438(8) 


0.470(5) 


0.499(6) 


0.501(4) 


a\ 


0.493(7) 


0.562(10) 


0.617(14) 


0.641(9) 


N 


0.690(21) 


0.742(9) 


0.773(10) 


0.781(5) 


N' 


0.731(16) 


0.802(15) 


0.872(13) 


0.912(16) 



tance rate with the new code we had to use a 
HMC time step size that was nearly half as large 
as that used previously. This is most likely due 
to some failure in our earlier program to adhere 
to the conventions of || . We are trying to resolve 
this discrepancy. 

In figure [j] we plot the staggered hadron masses 
of p and nucleon and their parity partners as a 
function of the dynamical quark mass using re- 
sults from simulation I, II and IV. To confirm the 
near degeneracy of the parity partners in the zero 
quark mass limit we have done a 3-parameter fit 
to the 6 data points for both parity partner pairs, 
forcing a common intercept in the chiral limit re- 
spectively. The fits have y^/dof = 1.2 for the 
N and N' data and x 2 /c?o/ = 0.7 for the p and 
a±, which strongly support our parity doubling 
observation. 

The degeneracy of parity partners implies a 
restoration of chiral symmetry breaking. Figure ^ 
shows {i^ip) as a function of the dynamical quark 
mass for both Nf — 4 and Nf = 2, where the 
four-flavor data (solid circle) is taken from run I, 
II and IV and the two-flavor data (open squares) 
from ||. A linear fit to the chiral condensate for 
four flavors (solid line) yields @ip) = 0.00278(33) 
in the zero-mass limit, which is more than three 



times smaller than (tpip) = 0.00854(17) for two 
flavors. This indicates significant weakening of 
chiral symmetry breaking as the number of fla- 
vors increases, although the four-flavor value is 
still about eight standard deviations from zero. 
The non-vanishing chiral condensate in the zero- 
mass limit suggests that the hadron mass split- 
tings between the parity partners will not be ex- 
actly zero although they are invisible in figure [j] 
because of the relatively large statistical errors on 
the hadron masses. 

In calculating the Wilson valence hadron prop- 
agators, we employ the operator i/'7573' ! /' f° r the 
a\ and '073V' f° r the p. Wall sources and Coulomb 
gauge are also used. In figure [3] the Wilson va- 
lence masses for the p and a\ measured in sim- 
ulation II, are plotted as a function of the va- 
lence hopping parameter k. Linear fits are per- 
formed to both sets of data. The n c is determined 
through the PCAC relation m 2 = C(1/k — 1/k c ). 
We can see that the parity partners split less and 
less as k approaches k c . Similar behavior is also 
observed in partially quenched staggered hadron 
masses from this and earlier runs [W . 

In conclusion, the new QCDSP simulations 
have reproduced our previous runs successfully. 
The new simulation with Nf = 4, j3 = 5.4 and 




Figure 1. Staggered masses for the p and N and 
their parity partners vs. dynamical quark mass 
for Nf = 4. The fits are three-parameter linear 
fits to six data points with common intercept. 



m = 0.015 on QCDSP produces results consistent 
with our previous potential discovery that chiral 
symmetry breaking is reduced substantially as Nf 
increases from 2 to 4. The Wilson valence hadron 
spectrum results are consistent with the hypoth- 
esis that this result is independent of lattice for- 
malism, although more simulations are needed. 
The possible effects of small volume remain an 
important concern when interpreting our results. 
To resolve this issue we need to simulate on lat- 
tices of larger spatial volumes. In fact, larger, 
24 3 x 32, full QCD simulations are now under- 
way. 
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Figure 2. (ipip) vs - dynamical quark mass. Filled 
circles are for four flavors and open squares are 
for two flavors. The solid line is a linear fit to the 
3 four flavor data points and the dashed line is a 
fit to the 4 two flavor data points. 




Ql , nip and vs. 1/k in run II 



Figure 3. m 
with Nf = 4 and md yn = 0.015. The filled circles 
are masses of p, open circles are masses of a± and 
diamonds are masses squared of n. The dashed 
vertical line shows the position of k c . 



